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(T) ' Abstract. We give algorithms for computing the singular moduli of suitable nonholomor- 

phic modular functions F(z). By combining the theory of isogeny volcanoes with a beauti- 
ful observation of Masser concerning the nonholomorphic Eisenstein series E%(z), we obtain 
CRT-based algorithms that compute the class polynomials Hj}(F;x), whose roots are the 
discriminant D singular moduli for F(z). By applying these results to a specific weak Maass 
form F p (z), we obtain a CRT-based algorithm for computing partition class polynomials, a 
sequence of polynomials whose traces give the partition numbers p(n). Under the GRH, the 
. expected running time of this algorithm is 0(n 5 / 2+0 W). Key to these results is a fast CRT- 

based algorithm for computing the classical modular polynomial $,„ (X, Y) that we obtain by 
extending the isogeny volcano approach previously developed for prime values of m. 



& \ 1. Introduction and Statement of results 



As usual, we let 



(l + 240£~ =1 £ d|n dV 



l ; 4 : - - .: -» < I I I ' I . 

CM ■ 1.1 j(z) := ± n 4 J = q- 1 + 744 + 196884g + 21493760g 2 + 

vo ; . . 

lO . be Klein's classical elliptic modular function on SL 2 (Z) (g := e 2mz throughout). The values 
of j(z) at imaginary quadratic arguments in the upper-half of the complex plane are known 
as singular moduli. Two examples are the Galois conjugates 

fl + y/^15\ -191025 - 85995^ , (\ + V^15\ -191025 + 85995\/5 
J { 2 )= 2 and J {-^—) = 2 ' 

These numbers play an important role in algebraic number theory. Indeed, they are algebraic 
integers that generate ring class field extensions of imaginary quadratic fields, and they are 
the j-invariants of elliptic curves with complex multiplication [3, 9, 10, 11]. 

The problem of computing singular moduli has a long history that dates back to the works 
of Kronecker, and is highlighted by famous calculations by Berwick [2] and Weber [35]. His- 
torically, these numbers have been difficult to compute. More recently, Gross and Zagier [21] 
determined the prime factorization of the absolute norm of suitable differences of singular 
moduli (further work in this direction has been carried out by Dorman [12, 13]), and Zagier 
[36] identified the algebraic traces of singular moduli as coefficients of half-integral weight 
modular forms. 
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Here we consider the problem of computing the minimal polynomials of singular moduli, 
the so-called Hilbert class polynomials. This problem has been the subject of much recent 
study. For example, Belding, Broker, Enge, and Lauter [1], and the third author [33], have 
provided efficient methods of computation that are based on the theory of elliptic curves with 
complex multiplication (CM). The basic approach in that work is simple. One uses theoretical 
facts about elliptic curves with CM to quickly compute the reductions of these polynomials 
modulo a set of suitable primes p, and one then compiles these reductions via the Chinese 
Remainder Theorem (CRT) to obtain the exact polynomials. The primes p are chosen in a 
way that facilitates the computation, and in particular, they split completely in the ring class 
field Kq of the imaginary quadratic order O associated to the singular moduli whose minimal 
polynomial one wishes to compute. Under the Generalized Riemann Hypothesis (GRH), this 
algorithm computes the discriminant D Hilbert class polynomial with an expected running 
time of 0(| D \)} In [6], the third author, together with Broker and Lauter, further developed 
these ideas to compute modular polynomials $^ using the theory of isogeny volcanoes. 

We extend these results to the setting of nonholomorphic modular functions such as 

(1-2) * W mZ) - WU 



QE 6 (z)j(z) 2y J §j{z){j{z)-l72%y 
where 

V 1 n=l d\n 

is the weight 2 nonholomorphic modular Eisenstein series, and where 

oo oo 

E A (z) := 1 + 240 J^l^ d V and E 6 {z) := 1 - 504 ^ 

n=l d\n n=l d\n 

are the usual weight 4 and weight 6 modular Eisenstein series. 

Remark. The function "f(z) plays an important role in this paper. We shall make use of an 
observation of Masser [27] which gives a description of its singular moduli in terms of the 
coefficients of certain representations of the classical modular polynomials. 

We first recall the setting of Heegner points on modular curves (see [20]). Let N > 1, 
and let D < be a quadratic discriminant coprime to N. The group T (N) acts on the 
discriminant D positive definite integral binary quadratic forms 

Q(X, Y) = [a, 6, c] := aX 2 + bXY + cY 2 

with N I a. This action preserves h (mod 2JV). Therefore, if /3 2 = D (mod 4JV), then it 
is natural to consider Qn,d,/3, the set of those discriminant D forms Q = [a,b,c] for which 
< a = (mod N) and b = (3 (mod 2N), and we may also consider the subset Qjv™ a 
obtained by restricting to primitive forms. The number of T (N) equivalence classes in Qn,d,/3 
is the Hurwitz-Kronecker class number H(D), and the natural map defines a bijection 

Qn,d,p/T (N) — > Q D /SL 2 (Z), 



1 We use the "soft" asymptotic notation 0{n) to denote bounds of the form 0(n log n). 
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where Qn is the set of discriminant D positive definite integral binary quadratic forms (see the 
proposition on p. 505 of [20]). This bijection also holds when restricting to primitive forms, 
in which case the number of r (iV) equivalence classes in Q^™^, and the number of SL 2 (Z) 

equivalence classes in Q^ im , is given by the class number h(D). 

For modular functions F(z) on SL^Z), our goal is to calculate the class polynomial 

(1.4) H D (F;x):= J] (x-F(a Q )), 

QeQP rim /SL 2 (Z) 

where ciq G H is a root of Q(x, 1) = 0. For modular functions F(z) on T (N) and discriminants 
D < coprime to N, our goal is to calculate the class polynomial 

(1.5) H D , p {F;x):= ]J {x-F(a Q )). 

We first consider the special case of the SL 2 (Z) nonholomorphic modular function 7(2), 
defined by (1.2), and we compute the Q-rational polynomials Hp^^x). 

Theorem 1.1. For discriminants D < —4 that are not of the form D = —3d 2 , Algorithm 1 
(see §3.1) computes ifo(7;x). Under the GRH, its expected running time is 0(\D\ 7 ^ 2 ) and it 
uses 0(\D\ 2 ) space. 

Remark. The discussion on p. 118 of [27], makes it clear how to modify Algorithm 1 to handle 
discriminants of the form D = —3d 2 . We expect that the bound on its expected running time 
can be improved to 0(\D\ 5 ^ 2 ) using tighter bounds on the size of the coefficients of Hd(^] x). 

A key building block of Algorithm 1 is a new algorithm to compute the classical modular 
polynomial $ m (X, Y), which parameterizes pairs of elliptic curves related by a cyclic isogeny 
of degree m. Here we extend the iosgeny volcano approach that was introduced in [6] to 
compute $ m for prime m so that we can now efficiently handle all values of m. The result is 
Algorithm 1.1 (see §3.1), which, under the GRH, computes $ m in 0(m 3 ) time. For suitable 
primes p it can compute $ m modulo p in 0(m 2 ) time (and space), which is crucial to the 
efficient implementation of Algorithm 1. 

Example. We have used Algorithm 1 to compute if^)(7,x) for D > —20000. Some small 
examples are listed below. 
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We extend Algorithm 1 to compute class polynomials for a large class of nonholomorphic 
modular functions. This class includes, for example, the SL 2 (Z)-function 

K( s _ 9RR E*(z)E,(z)E 6 (z) + ?>E A {zf + 2E 6 (z) 
( j E 4 (z)*-E 6 (zy 

considered by Zagier (see §9 of [36]) in his famous paper on traces of singular moduli. More 
generally, it includes suitable modular functions F(z) of the form 

F(z) :=9_ a o9_4...^_ 2fc (y), 

where $(z) is a weight 2 — 2k weakly holomorphic modular form on r (iV) whose Fourier 
expansions at cusps are algebraic. Here the differential operator dh, which maps weight h 
modular forms to weight h + 2 modular forms, is defined by 

(1.6) d h ■ 



2ni dz 47rlm(2;) 

We consider a specific class of such modular functions. Let O be the imaginary quadratic 
order with discriminant D, ring class field Kq, and fraction field K = Q(yD). Let c\ and c 2 
denote fixed positive integers. Let F(z) be a modular function (i.e. weight 0) for To(N) that 
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can be written in the form 



where each A n 6 Q(j) is a rational function of j(z). The function F(z) is said to be good for 
a discriminant D < coprime to N if it satisfies the following: 

(1) Each A n (a Q ) lies in K Q for all Q G Q^ rim /SL 2 (Z). 

(2) The polynomial Ci\D\ C2h H£>(F; x) has integer coefficients. 

Remark. The class of good modular functions includes many nonholomorphic modular func- 
tions, such as 7(2), and it includes meromorphic modular functions which may have poles in 
the upper half plane (poles at CM points are excluded by condition (1)). 

For good modular functions, we obtain the following general result. 

Theorem 1.2. For discriminants D < —4 not of the form D = —3d 2 , Algorithm 2 (see §3.2) 
computes class polynomials for good modular functions F(z). Under the GRH, its expected 
running time is 0{\D\^ 2 ), and it uses 0(\D\ 2 ) space. 

In fact, Algorithm 2 can be readily adapted to treat modular functions of the form F(z) = 
A n {z) r y{z) n where the coefficient functions A n (z) do not necessarily lie in Q(j), using the 
techniques developed by Enge and the third author in [16]. As an example, we apply Algo- 
rithm 2 to obtain a CRT-based algorithm for computing the "partition polynomials" defined 
by the first two authors in [8]. These are essentially the class polynomials of the r (6) non- 
holomorphic modular function 

(1.7) F M := -<MP W , = (l - ^L_) r > + + (29 + ^|_) , + 
where P{z) is the weight —2 weakly holomorphic modular form 

;= 1 . aw - ay - m + 6g, (fe) = , _ 10 _ 29q + 

2 i]{z) 2 r]{2z) 2 r](3z) 2 r]{6z) 2 
These polynomials are defined as 

(1.8) Hr\x):= J] (x-P(a Q )). 

Q£S6,l-24n,l 

In contrast with (1.4) and (1.5), we stress that the roots of these polynomials include singular 
moduli for imprimitive forms (if any). The interest in these polynomials arises from the fact 
that (see Theorem 1.1 of [8]) 

(1.9) Hl w \x) = x h{1 ~ 24n) - (24n - l^n)^ 1 - 24 ^ 1 + . . . , 

where p(n) is the usual partition function. Since the roots P(cxq) are algebraic numbers that 
lie in the usual discriminant 1 — 24n ring class field, we have the following finite algebraic 
formula 

p( n ) = V P(a ). 

yy 1 24n-l ^ V Q! 

Q£Q6,l-24n,l 

Remark. By the work of the first two authors [8], combined with recent results by Larson and 
Rolen [26], it is known that each (24n — l)P(ag) is an algebraic integer. 
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As a consequence of Theorem 1.2, we obtain the following result. 

Theorem 1.3. For all positive integers n, Algorithm 3 (see §3.3) computes H^ vt {x). Under 
the GRH, its expected running time is 0(n 5//2 ) and it uses 0(n 2 ) space. 

Remark. For the simpler task of computing individual values of p(n), an efficient implemen- 
tation of Rademacher's formula such as the one given in [23] is both asymptotically and 
practically faster than Algorithm 3. 

Example. We have used Algorithm 3 to compute H^ &vt (x) for n < 750. Some small examples 
are listed below. 



n (24n - l)p(n) H^ rt (x) 

1 23 x 3 - 23a; 2 + ^ x - 419 

2 94 x 5 _g 4x 4 + ™ x 3_ 65838x 2 + 10™i76 x+ 14^)23 

3 213 x 7 - 213a; 6 + ™^x 5 - 723721a; 4 + 44648 ^f 2886 x 3 

I 9188934683 2 _i_ 166629520876208 , 2791651635293 
< 71 X ~r yj3 X T" 

4 475 x s - 475a; 7 + - 9455070a; 5 + 39495 ^f 9743 x 4 

95 95 z 

_ 97215753021 3 i 9776785708507683 2 
'X T" oc;3 



19 1 95 3 



53144327916296_ x _ 134884469547631 



19 2 5 4 -19 



This paper is organized as follows. In §2 we recall essential facts about elliptic curves with 
complex multiplication, and singular moduli for modular forms and certain nonholomorphic 
modular functions. In §3 we use these results to derive our algorithms. In §4 we conclude 
with a detailed example of the execution of Algorithm 3 for n = 1 and n = 24. 

2. NUTS AND BOLTS 

We begin with some preliminaries on elliptic curves with complex multiplication and singular 
moduli for suitable modular functions. 

2.1. Elliptic curves with complex multiplication. We recall some standard facts from 
the theory of complex multiplication, referring to [9, 25, 30] for proofs and further background. 
Let O be an imaginary quadratic order, identified by its discriminant D. The j-invariant of the 
lattice O is an algebraic integer whose minimal polynomial is the Hilbert class polynomial Ho- 
If a is an invertible C-ideal (including a = O), then the torus C/a corresponds to an elliptic 
curve E/C with complex multiplication (CM) by O, meaning that its endomorphism ring 
End(.E) is isomorphic to O, and every such curve arises in this fashion. Equivalent ideals 
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yield isomorphic elliptic curves, and this gives a bijection between the ideal class group cl(0) 
and the set 

(2.1) EUo(C) := {j(E/C) : End(£) = O}, 

the j-invariants of the elliptic curves defined over C with CM by O. We then have 

(2.2) H D (x):= l[ (x- Jt )=H D (j;x). 

j<6EUo(C) 

The splitting field of Hp over K = Q(y/~D) is the ring class field Kq. It is an abelian extension 
whose Galois group is isomorphic to cl(0), via the Artin map. 

This isomorphism can be made explicit via isogenies. Let E/C be an elliptic curve with CM 
by O and let a be an invertible C-ideal. There is a uniquely determined separable isogeny 
whose kernel is the subgroup of points annihilated by every endomorphism in a C O 
End(-E). The image of this isogeny is an elliptic curve that also has CM by O, and this 
defines an action of the ideal group of O on the set Ello(C). Principal ideals act trivially, and 
the induced action of the class group is regular. Thus Ello(C) is a principal homogeneous 
space, a torsor, for the finite abelian group cl(0). 

If p is a (rational) prime that splits completely in Kq, equivalently, for p > 3, a prime 
satisfying the norm equation 

(2.3) Ap = t 2 - v 2 D 

for some nonzero integers t and v, then Hp splits completely in ¥ p [x] and its roots form the 
set 

(2.4) EU (F P ) := {j(E/¥ p ) : End(£) = O}, 

Conversely, every ordinary (not supersingular) elliptic curve E/¥ p has CM by some imaginary 
quadratic order O in which the Frobenius endomorphism corresponds to an element of norm p 
and trace t. 

2.2. Modular polynomials via isogeny volcanoes. For each positive integer m, the clas- 
sical modular polynomial $ m is the minimal polynomial of the function j(mz) over the field 
C(j). As a polynomial in two variables, $ m G Z[X, Y] is symmetric in X and Y. If E/k is 
an elliptic curve and N is prime to the characteristic of k, then the roots of <& m (j(E), Y) are 
precisely the j-invariants of the elliptic curves that are related to E by a cyclic m-isogeny; 
see [25] for these and other properties of $ m . 

For distinct primes i and p, we define the graph of i-isogenies Ge(¥ p ), with vertex set ¥ p 
and edges (ji, J2) present if and only if $^(ji,j2) = 0. Ignoring the connected components 
of and 1728, the ordinary components of Ge(¥ p ) are l-volcanoes [17, 24], a term we take 
to include cycles as a special case; see [34] for further details on isogeny volcanoes. In this 
paper we focus on £- volcanoes of a special form, for which we can compute $^ mod p in a 
particularly efficient way, using [6, Alg. 2.1]. 

Let O be an order in an imaginary quadratic field K with maximal order Ok, and let i be 
an odd prime not dividing [Ok '■ O}. Assume D = disc(C) < —4. Suppose p is a prime of the 
form 4p = t 2 — £ 2 v 2 D with p = 1 mod t and £ J( v; equivalently, p splits completely in the ray 
class field of conductor £ for O and does not split completely in the ring class field of the order 
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with index £ 2 in O. Then the components of Ge(¥ p ) that intersect Ello(F p ) are isomorphic 
£- volcanoes with two levels: the surface, whose vertices lie in Elle>(Fp), and the floor, whose 
vertices lie in Elle>/(F P ), where O' is the order of index £ in O. Each vertex on the surface is 
connected to 1 + (-j) = 0, 1 or 2 siblings on the surface, and £ — (y) children on the floor. An 
example with £ = 7 is shown below: 




Provided h(0) > £ + 2, this set of ^-volcanoes contains enough information to completely 
determine <3>^ mod p. This is the basis of the algorithm in [6, Alg. 2.1] to compute modp, 
which we make use of here. Selecting a sufficiently large set of such primes p allows one to 
compute over Z (via the CRT), or modulo an arbitrary integer M (via the explicit CRT). 
Our requirements for the order O and the primes p are summarized in the definition below. 

Definition 2.1. Let i > 2 be prime, and let c > 1 be an absolute constant independent of £. 
An imaginary quadratic order O is said to be suitable for I if I / [Ok '■ O] and £ + 2 < h(0) < 
c£. A prime p is then said to be suitable for £ and O if p = 1 mod £ and Ap = t 2 — £ 2 v 2 disc(O), 
for some t, v G Z with £ J( v. 

The definition of suitability above is weaker than that used in [6], but this only impacts 
logarithmic factors in the running time that are hidden by our soft asymptotic notation. 

2.3. Selecting primes with the GRH. In order to apply the isogeny volcano method to 
compute (or the polynomials H^(F;x) we wish to compute), we need a sufficiently large 
set S of suitable primes p. We deem S to be sufficiently large whenever ^ pg5 logp > B + log 2, 
where B is an upper bound on the logarithmic height of the integer coefficients that we wish to 
compute with the CRT. 2 For &i{X, Y) = a ijX l Y^ , we may bound ht(<£v) := logmaxjj |ay| 
using 

(2.5) ht(^) < Q£\og£+ 1S£, 

as proved in [7]. 

Heuristically (and in practice), it is easy to construct the set S. Given an order O of 
discriminant D suitable for £, we fix v = 2 if D = 1 mod 8 and v — 1 otherwise, and for 
increasing t = 2 mod £ of correct parity we test whether p = (t 2 — v 2 £ 2 D)/4 is prime. We add 
each such prime to S, and stop when S is sufficiently large. 

Unfortunately, we cannot prove that this method will find any primes, even under the GRH. 
Instead, we use Algorithm 6.2 in [6], which picks an upper bound x and generates random 
integers t and v in suitable intervals to obtain candidate primes p = (t 2 — v 2 £ 2 D)/A < x that 
are then tested for primality. The algorithm periodically increases x, so its expected running 
time is 0(B 1+e ), even without the GRH. 

Under the GRH, there are effective constants C\, c 2 > such that x > C\£ 6 log 4 £ guarantees 
at least c 2 £ 3 \og 3 £ suitable primes less than x, by [6, Thm. 4.4]. Asymptotically, this is far 
more than the 0(£) primes we need to compute $£. We note that S contains 0(B j log B) 
primes (unconditionally), and under the GRH we have logp = 0(log-B + log/) for all p G S. 



'When the coefficients are rational numbers that are not integers, we first clear denominators. 
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2.4. Modular singular moduli. The results from the previous section can be cast in terms 
of the CM values of the j-function. Indeed, we have the following classical theorem (for 
example, see [3, 9]) which summarizes some of the most important properties of singular 
moduli for Klein's j-function. 

Theorem 2.2. Suppose that Q = ax 2 + bxy + cy 2 is a primitive positive definite binary 
quadratic form with discriminant D = b 2 — 4ac < ; and let cvq G H be the point for which 
Q{c(q, 1) = 0. Then the following are true: 

(1) The singular modulus j(q;q) is an algebraic integer whose minimal polynomial has 
degree equal to the class number h(D). 

(2) The Galois orbit of j (cxq) consists of the j(z)-singular moduli associated to the h(D) 
equivalence classes in Q^ im /SL 2 (Z). 

(3) If K = Q{\/~D), then the discriminant D singular moduli are conjugate over K . More- 
over, K(j(a.Q)) is the ring class field of the quadratic order of discriminant D; in the 
case that D is a fundamental discriminant, K{j{ctQ)) is the Hilbert class field of K. 

Theorem 2.2 and the properties of the weight 2 nonholomorphic Eisenstein series E%{z) 
at CM points shall play a central role in the construction of the algorithms described in the 
next section. To this end, we make use of the special nonholomorphic function j(z) defined 



Masser nicely observed that the singular moduli for 7(2) can be computed using the singular 
moduli for j(z) and certain expressions for modular polynomials. Here we make this precise 
for discriminants D < —4 that are not of the form D = —3d 2 . 

Remark. Masser explains how to handle discriminants D = —3d 2 ; see p. 118 of [27]. 

To state his observation, we let O be the imaginary quadratic order of discriminant D, and 
let {Qi, . . . , Qh} be a set of representatives for Q^ im /SL 2 (Z) ~ c\(0), where h = h(D) is the 
class number. To simplify notation, we use $£> to denote the classical modular polynomial 
$|£)|(X, Y). For any Q = Qi we may write in the form 



which satisfies if} {mi) = 0{m log log m); see [32]. The coefficients (3^ u = (3^ u {olq) are algebraic 
integers that lie in the ring class field Kq, and we have (3^ = (3 Uifl (by the symmetry of $£>). 
Masser [27, p. 118] gives the following formula for 7(ckq). 

Lemma 2.3. Assuming the notation and hypotheses above, we have 



in (1.2). 



(2.6) 




0<fi,v<n 



where n — ip{\D\) is determined by the Dedekind ^-function 




p\m 



(2.8) 



7(«q) 



2/3 , 2 («q) - /3i,i(o!q) 
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Two remarks. 

1) Masser proves (see [27, Lem. A2]) that /3o,i(c*q) is nonzero. 

2) From (2.6), one finds that 

M<x Q ) = [Y]$ D (j{a Q ),Y +j(a Q )), 
(2-9) A,i(«q) = [Y]& D (j(a Q ) 7 Y + j(a Q )) 7 

Molq) = [Y 2 ]^ D (j(a Q ),Y + j(a Q )), 

where & D (X,Y) = -^$ D (X,Y), and for any polynomial f(Y), the notation \Y k ]f(Y) indi- 
cates the coefficient of Y k in f(Y). 

3. The Algorithms 

Here we apply and extend the results in §2 to derive our algorithms. 

3.1. Algorithm 1. We now give an algorithm to compute the class polynomial Hd(j;x), 
where j(z) is the nonholomorphic modular function defined in (1.2) and D is an imaginary 
quadratic discriminant. In order to simplify the exposition as above, we shall assume D < —4 
and that D is not of the form D = —3d 2 ; these special discriminants are in principle no more 
difficult to handle than the general case, but the details are more involved; see [27, p. 118]. 3 
To make use of Lemma 2.3, we need to compute the singular moduli j(ag). These shall be 
obtained as the roots of the Hilbert class polynomial Hd(x). Thus if we know $d and Hp, 
then we can apply Lemma 2.3 to compute 

H D (-f;x)= Yl ( x -l( a Q))- 

QeQP rim /SL 2 (Z) 

Using algorithms for fast multipoint polynomial evaluation and fast integer arithmetic (see 
[18], for example), this yields an algorithm that computes H D (^;x) in 0(\D\ 3 ) expected time 
using 0(|D| 3 ) space, under the GRH. However, this approach is quite memory intensive and 
quickly becomes impractical, even for moderate values of D. As an alternative, we give a 
CRT-based algorithm that uses 0(\D\ 7 / 2 ) expected time and 0(\D\ 2 ) space, under the GRH. 

It is clear from equations (2.6) and (2.8) that the coefficients of Hd(j',x) lie in Q; indeed, 
the coefficients of $^ are integers, as are the elementary symmetric functions of the j(«q), 
which are the coefficients of the Hilbert class polynomial Hd(x). If we let 

(3.1) S:= J] Po,i{<*q), 

QeQP rim /SL 2 (Z) 

then 5 G Z is divisible by the denominator of every coefficient of Ho{'j',x) and 5H£ ) {'y\x) G 
7j[x]. We now present the algorithm. 

Algorithm 1 

Input: An imaginary quadratic discriminant D that is not special. 

3 We note that the discriminants D = 1 — 24n needed to compute H^ alt (x) are not special. 
4 As remarked in the introduction, we expect this running time can be improved, possibly to 0(\D\ 5 / 2 ), by 
obtaining tighter bounds on the coefficients of iJo(7; x). 
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Output: The polynomial Hd(j; x) G Q[x\. 

1. Pick an order O suitable for |D|, and a set 5* of primes suitable for \D\ and O (see Def. 3.2), 
using the bound B 1 (D) given in (3.3) below. 

2. Compute the Hilbert class polynomial Hp G Z[x] using [33, Alg. 2]. 

3. For each prime p G S: 

a. Compute $ = $d mod p using Algorithm 1.1 below. 

b. Compute $' = K) mod p. 

c. Compute the roots ji, . . . ,jh G ¥ p of ifp mod p. 

d. Compute fc (F) = $(j k ,Y) and </>' fe (Y) = &(Jk,Y) for a11 h using [18, Alg. 10.7]. 

e. For each j k , compute /0i,i> an d A),2 using fe and via (2.9), 
and then compute ^ k = (2/3 ,2 — / /?o,i- 

f. Compute 5 = Yl k /3 ,i and f(x) = 5 Yl k {x - j k ). 

g. Save f(x) mod p and 5 mod p. 

4. Use the CRT to recover f(x) = 5Hd(j', x) G Z[x] and 5 G Z. 

5. Output tf D ( 7 ;a;) = G Q[x]. 

Let B$(D) denote an upper bound on ht ($d); when |D| is prime we may use the bound 
Bq>(D) = 6\D\ + 18\D\ log \D\ from (2.5), and otherwise we may derive such a bound by 
expressing $^ in terms of modular polynomials of prime level, as in [9, Thm. 13.14]. We use 

(3.2) M(D) := log(exp(7rv^Dj) + 2114.567) 

to bound log |i(«g)|, for any Q G Q^ im ; see [14, p. 1094], for example. 
We now define 

(3.3) B 7 (D) := (h(D) + l)(41ogty(\D\) + l)+2i;(\D\)M(D) + Bg,(D) + 2\ 

Lemma 3.1. Let D be an imaginary quadratic discriminant that is not special. Then we have 
the bounds 5 < B^D) and ht(5H D ( T , x)) < B y (D), and 5 = 6(|Z>| 3 / 2 ). 

Proof. For any Q G Q^ im , each coefficient c of the univariate polynomial <f>(Y) = $dU(o!q), Y) 
is a polynomial of degree i/;(\D\) m 3( a o) with coefficients bounded by ht($£>) < B$(D). It 
follows that 

(3.4) ht(0) < ip(\D\)M(D) + B*(D) + log(^(|D|) + 1). 

From (2.9), we know that A),i(o;q) is the linear coefficient of <fi(Y + j(oiq)), where 4>(Y) has 
degree i[)(\D\), and this implies 

(3.5) log \Po,i(a Q )\ < 2 bg(V(|D| + 1) + ht(0) + ijj{\D\)M 

for all Q G Qjp im . Substituting (3.4) into 3.5 and applying the bound to each of the h(D) 
factors in the product 5 = Yin A>,i( a Q) yields 

(3.6) log \8\ < h(D) (3 log(V(|D|) + 1) + 2iJ>(\D\)M(D) + B$(D)) < B 7 (D). 
A calculation completely analogous to that used in (3.5) yields 

(3.7) log |2/3 , 2 (a Q ) + /?i,i(«q)| < 3 log(^(|D| + 1) + ht(0) +^(\D\)M + 2, 



12 



JAN HENDRIK BRUINIER, KEN ONO, AND ANDREW V. SUTHERLAND 



for all Q G Q^ im . The absolute values of the numerators of the coefficients of Hd(Ti x )> 
which has degree h(D), have logarithms that exceed the bound in (3.7) by at most \og2 h ( D \ 
Combining this with (3.6) yields 

(3.8) H(6H D (j; x)) < log \5\ + 3 bg(^(|D| + 1) + ht(0) + ip(\D\)M + +h(D) log2 + 2, 

and it is then easy to check that plugging (3.4) and (3.6) into the RHS yields an expression 
that is bounded by B^(D) as defined in (3.3). The asymptotic bound B^(D) = 0(\D\ 3 / 2 ) 
follows immediately from the bounds h(D) = 0{\D\ 1 / 2 ), i/)(\D\) = 0(\D\), M(D) = 0(\D\V 2 ) 
and B$(D) = 0(\D\). □ 

For an odd prime £, given a suitable order O and a suitable prime p, the isogeny volcano 
algorithm of [6, Alg. 2.1] computes $^ in 0(£ 2 ) time, provided that logp = 0{\og£). Here we 
extend this result to any integer m > 1. We first note that 

$ 2 (X, Y) =X 3 + Y 3 - X 2 Y 2 + 1488(X 2 y + XY 2 ) - 162000(X 2 + Y 2 ) 
+ 40773375XF + 8748000000 (X + Y) - 157464000000000, 

and extend Definition 2.1 to composite integers m. 

Definition 3.2. Let m > 1 be an integer, let £ be the largest prime divisor of m. An imaginary 
quadratic order O is said to be suitable for m if if)(m) + 1 < h(0) < 3ip(m) and £ / [Ok '■ O]. 
A prime p is suitable for m and O if p = 1 mod £ and 4p = t 2 — £ 2 v 2 disc(O) for some t, v G Z 
with £ I v. 

Algorithm 1.1 

Input: An integer m > 1, an order O suitable for m, and a prime p suitable for m and O. 
Output: The modular polynomial $ m mod p. 

1. If m = 2 then output $2 mod p and terminate. 

2. If m is an odd prime then compute $ m mod p via [6, Alg. 2.1] and terminate. 

3. Compute $£ for each prime £ < ^/m dividing m. 

4. Compute the Hilbert class polynomial H D , where D = disc(C), via [33, Alg 2]. 

5. Compute the roots ji,...,j h G ¥ p of H D mod p and let S = {ji, ■ ■ ■ , j n }, where n = 
ip(m) + 1. 

6. If m has a prime divisor £ > \fm then compute n sets of j'-invariants S®, . . . , 5°, where 
5° = {J G F p : $£ (jj, J) = 0}, using the isogeny volcano method. 

Otherwise, let £0 = 1, and let S® = {ji} for 1 < i < n. 

7. Let £\ < £ 2 < ■ ■ ■ < l T be the primes whose product is m/£o. 
For 1 < d < r do the following: 

a. If £d 7^ then compute Sf = {j G ¥ p : $^(jfc, J) for some ^ for 1 < i < n. 

b. If £ d = £ d . x then compute Sf = {j G ¥ p : $ £i ( Jfc , j) for some Jfc G Sf-^Sf' 2 } for 1 < 
i <n. 

8. For 1 < z < n compute <fii(X) = $ m (X, jj) = Ylt=o a ikX k as the product rijes"-(^ — j)- 

9. For < k < n interpolate the polynomial fk{X) of degree less than n for which fk{ji) = &ik- 
10. Output $ m (X,Y) = YlS fkY k mod p. 
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Note that step 6 does not use $£ to compute S®, it uses the isogeny volcano method detailed 
in [6, §6], whereas step 7 uses the (smaller) polynomials $£ computed in step 3. When 
computing $ m mod p for many primes p (as in Algorithm 1), the polynomials Hd G 1>[x\ and 
G Z[X, Y] computed in steps 3 and 4 may be computed just once and reused, since they 
do not depend on p. 

Lemma 3.3. Algorithm 1.1 correctly computes $ m (X, Y) mod p. Under the GRH, its expected 
running time is 0(m 2 ), provided that logp = O(logm). 

Proof. For each of the j-invariants ji G S, the set Sf contains all the j-invariants J for which 
&m d (ji,j) = 0, where = nf=o^ d - This follows from the defining property of $ m (X, Y) (it 
parameterizes cyclic m-isogenies) and the fact that every (separable) cyclic isogeny can be 
expressed as a product of cyclic isogenics of prime degree (note that p is distinct from all the 
£i, since p = 1 mod £{). Thus we have (j>i{X) = $ m (X, jj) in step 8, and the n = ip(m) + 1 
distinct values of ji G S are sufficient to uniquely determine the coefficients of $ m mod p in 
step 9. 

We now bound the complexity of each step, assuming the GRH and that logp = O(logm). 
Step 1 takes 0(1) time and step 2 takes 0(m 2 ) time, by [6, Thm. 6.5]. Computing O(logm) 
modular polynomials $£ with £ < y/m in step 3 takes 0(m 3//2 ) expected time, by [6, Thm. 1]. 
We have h = h(D) = 0(m), since O is suitable for m, and since h(D) = 0(|D| 1//2 ), we have 
\D\ = 0{m 2 ) and this bounds the cost of step 3, by [33, Thm. 1]. Using standard probabilistic 
algorithms for root-finding, step 5 takes 0(h) expected time, which is 0(m). 

The cost of step 6 is bounded by O(h£ +£l), which is 0(m 2 ); this follows from the proof of [6, 
Thm. 6.5]. Step 7 consists of root-finding operations in F p whose expected complexity is softly- 
linear in the number of roots, ignoring factors of logp (every polynomial under consideration 
splits completely ¥ p [x] by virtue of the suitability of O and p). The total number of roots 
computed in step 7 is 0(ip(m) 2 ), hence the total cost is 0{m 2 ) expected time. 

Using standard algorithms for fast arithmetic and polynomial interpolation [18], the cost 
of steps 8 and 9 are both bounded by 0(ip(m) 2 ), which is 0(m 2 ). Thus every step has an 
expected running time bounded by 0(m 2 ). □ 

Corollary 3.4. Under the GRH, for any integer m > 1 the modular polynomial $ m can be 
computed in 0(m 3 ) expected time. 

Proof. An explicit 0{m) bound on the height of $ m can be derived from [9, Prop. 13.14] using 
the height bounds for for primes £\m given in (2.5). An order O suitable for m can be 
obtained from the family of suitable orders given in [6, Ex. 4.3] and a sufficiently large set S 
of primes p suitable for m and O can be selected using [6, Alg. 6.2]. Under the GRH, these 
primes satisfy logp = O(logm) and the corollary then follows from Lemma 3.3 and standard 
bounds on the time for fast Chinese remaindering [18, §10.3]. □ 

Remark. An algorithm to compute $ m using floating point approximations appears in [15] 
with a running time that is also 0(m 3 ), but the correctness of this algorithm and the bound 
on its running time both depend on a heuristic assumption regarding the precision needed to 
avoid rounding errors. We note that Algorithm 1.2 is faster in practice, its output is provably 
correct, and the bound on its expected running time depends only on the GRH. 
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Corollary 3.5. Let D be an imaginary quadratic discriminant that is not special, let O be 
an order suitable for \D\ with h(D) = 0(\D\), and let p be a prime suitable for \D\ and O. 
Under the GRH, the polynomial #0(7; x) mod p can be computed in 0(\D\ 2 ) expected time. 



We now prove Theorem 1.1 given in the introduction, which we restate here. 

Theorem 1.1. Algorithm 1 computes i/o(7; x). Under the GRH, its expected running time 
is 0(\D\ 7 / 2 ) and it uses 0(\D\ 2 ) space. 

Proof. The correctness of Algorithm 1 follows from the discussion preceding the algorithm, 
which shows how to compute Hui^y-x) in terms of <3>d and the Hilbert class polynomial 
Hp, the correctness of the algorithm used to compute Hp [33, Thm. 1], the correctness of 
Algorithm 1.1 used to compute $£> modp (Lemma 3.3), and the validity of the bound B 7 (D) 
on the logarithmic height of Hd(^;x) (Lemma 3.1). 

For the time and space bounds, we now assume the GRH. We first note that, as explained 
in §2.3, we can select the set of primes S in step 1 in 0(B 7 (D) 1+e ) time, which is 0(\D\ 3/2+e ) 
for any e > 0, since B 7 {D) is 0(\D\ 3 / 2 ), by Lemma 3.1, and the primes p £ S all satisfy 
logp = 0(log \ D\). The time to compute the Hilbert class polynomial Ho in step 2 is 0(\D\), 
by [33, Thm. 1], and its size is 0{\D\). The set S has cardinality 0{B y (D)), which is 0(\D\ 3 / 2 ), 
and each iteration of step 3 takes 0(|D| 2 ) expected time: this follows from Lemma 3.3, the 
time for fast multipoint polynomial evaluation [18, Cor. 10.8], and standard bounds on the 
complexity of fast arithmetic in F p [x]. Thus the total expected time for step 3 is 0(\D\ 7 / 2 ). 

The space used in each iteration of step 3 must be bounded by 0(\D\ 2 ), since this bounds 
the time, and the total size of the values f(x) mod p and 5modp saved is bounded by 
0{h{D)B 1 {D)) 1 which is 0(|£>| 2 ). Finally, with fast Chinese remaindering [18, Alg. 10.22], 
the cost of step 4 is softly-linear in the total size of the coefficients of <5-£/d(7; x) and 5, which 



3.2. Algorithm 2. We now give an algorithm to compute the class polynomial Hd(F; x) for 
a good modular function F(z) = ^2 A n (z) , ~f(z) n , as defined in the introduction. We assume 
that each A n (z) is written in the form A n (z) = r n (j(z)), where r n £ Z(x). 

Algorithm 2 

Input: An imaginary quadratic discriminant D that is not special. 
Output: The polynomial Hd(F;x) £ Q[x\. 

1. Pick an order O lying in the order of discriminant D that is also suitable for \D\, and a 
set S of primes suitable for \D\ and O (see Def. 3.2) such that no prime in S divides the 
denominator of any of the r n (x), using the height bound Bp(D) (discussed below). 

2. Compute the Hilbert class polynomial Hd £ "L\x\ using [33, Alg. 2]. 

3. For each prime p £ S: 

a. Compute modp using Algorithm 1.1. 

b. Compute the roots ji, . . . ,jh £ F p of Hd mod p. 

c. For each jk do the following: 



Proof. Apply Lemma 3.3 to step 3 of Algorithm 1. 



□ 



is 0(\D\ 2 ). 



□ 
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i. Compute 7^ = (2/3 ,2 — /3i,i)/Ad,i modp as in Algorithm 1. 

ii. Compute F k = Y, r n{jk)lk m odp. 

d. Compute f(x) = c\\D\ C2h Y\ k (x — F k ) mod p. 

e. Save f(x) mod p. 

4. Use the CRT to recover f(x) = Cl \D\ C2h ^ H D (F; x) E Z[x]. 

5. Output H D {F-x) = ±\D\- C2h ^f (x) e Q[x). 



The bound Bp{D) used in step 1 is an upper bound on ht {c x \D\ C2h ^ H D (F] x)), which 

the next result shows is 0(|D| 1//2 ). Explicit computation of Bp(D) depends on the particular 
functions A n (z); bounds on the heights of the class polynomials H£>(A n ;x) can be readily 
derived from the functions r n (x) and known bounds on the height of the Hilbert class poly- 
nomial Hp; see Lemma 8 in [33], for example. From these, one can derive an explicit bound 
Bp(D) on the height of Hd(F;x); see Lemma 3.9 in the next section for an example. In 
general, the following lemma gives us an asymptotic bound for Bp(D) that suffices to bound 
the complexity of Algorithm 2. 

Lemma 3.6. For all non-special imaginary quadratic discriminants D we have 

ht ( Cl \D\ C2h ^H D (F;x)) = 0(\D\ 1/2 ). 

Proof. The proof follows as in the proof of Lemma 8 of [33]. One only needs to take care of 
the dependence of the summand j^j- - in the definition of E%(z) which in turn appears in 

the definition of j(z). We leave these details to the reader. □ 

We now prove Theorem 1.2 given in the introduction, which we restate here. 

Theorem 1.2. For discriminants D < —4 not of the form D = —3d 2 , Algorithm 2 computes 
H£>(F;x) for good modular functions F(z). Under the GRH, its expected running time is 
d(\D\ 5 ' 2 ) and it uses 0(\D\ 2 ) space. 

Proof. The correctness of Algorithm 2 is clear. We now bound its complexity, under the 
GRH. By Lemma 3.6, the set S contains 0(|D| 1//2 ) primes, and as described in §2.3, we can 
construct S in 0(\D\ l l 2 ) expected time. The expected time to compute Hjj(X) in step 2 is 
0(\D\), by [33, Thm. 1]. Each of the primes p G S satisfies \ogp = 0(log \D\), and therefore 
the expected time for step 3a is 0(\D\ 2 ), by Lemma 3.3. This dominates the cost of steps 3b 
and 3c, and the total expected time for step 3 is thus 0(\D\ 5 ^ 2 ), which dominates the expected 
running time of the entire algorithm. The space complexity of step 3 is dominated by the size 
of $d mod p, which is 0(\D\ 2 ). □ 

3.3. Algorithm 3. We now give an algorithm to compute the partition polynomial 

HT\x) := J] (x- P{a Q )) 

Q£S6,l-24n,l 

defined in (1.8). We do this by expressing H^ art (x) as a product of class polynomials 

(3.9) H D , p {P;x):= J] (x-P(a Q )). 

Q6Sri m «/ro(6) 
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Lemma 3.7. For a positive integer n, let D = 1 — 24n. Then 

u 2 \D 

where the class polynomials on the right and side are defined by (3.9), and e(u) = 1 if u = 
±1 mod 12 and e(u) = — 1 otherwise. 

Proof. Using [20], p. 505 equation (1), we obtain 

(3-10) Hr\x)= J] tf (1 _ 24n)/u2Ai (P;x), 

« 2 |l-24n 

where j3 u G Z/12Z denotes the unique residue such that (3 U ■ u = 1 (mod 12), equivalently, 
(3 U = u (mod 12). 

Let A < be any discriminant such that A = 1 (mod 24). Since P is invariant under the 
Atkin-Lehner involution We, we have 

H A ^(P; x) = H AA (P;x). 

Since P is taken to its negative under the Atkin-Lehner involution W$, we have 

H A ,± 5 (P;x) = H A ,x(-P;x) = (-l) h ^H A>1 (P;-x). 

Putting this into (3.10), we obtain the assertion. □ 

For the remainder of this section (and also in Section 4), we shall abuse notation and drop 
the dependence on /3 for modular functions on r (6). In every case we will have (3 = 1. For ex- 
ample, we let Hd(P\ x) := Hjj t i(P; x) for convenience. We cannot directly apply Algorithm 2 
to compute Hd(P; x) because the function P(z) does not satisfy all of our requirements for a 
good modular function; some minor changes are required, as we now explain. 

As shown by Larson and Rolen, the function P(z) may be decomposed as 

(3.H) P(z)=A(z) + B(z) 7 (z), 

where A(z) = A(z)j(z)(j(z) — 1728) and B(z) are weakly holomorphic modular functions for 
To(6); see Lemma 2.2 in [26]. The expression for P(z) in (3.11) does not satisfy our definition 
of a good modular function F(z) because A(z) and B(z) are not rational functions of j(z). 
However our two key requirements are satisfied: for discriminants D of the form 1 — 24n, 
the values of A(z) and B(z) at CM points lie in the ring class field and the polynomial 
\D\H D (P;x) has integer coefficients (so c\ = c 2 = 1). 
For each of the functions g = A,B, Larson and Rolen compute explicit polynomials 

* g (X;z) :=l[(X-g(az)), 

a 

where the product ranges over right coset representatives a for To (6) in SL2(Z). The polyno- 
mials ty g may be expressed as polynomials in X whose coefficients are integer polynomials in 
j(z), and we regard them as elements of Z[A, J]; see Appendix A of [26] for the exact values 
of $ 2 an d *&b- Here each occurrence of j(z) is replaced by the indeterminate J. While A(z) 
and B(z) are not rational functions of j(z), we note that the curves defined by ^^(X, J) and 
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^b(X, J) both have genus (and thus admit a rational parametrization, although we shall 
not make explicit use of this fact). 

It follows that, at least for discriminants prime to 6, the CM values of A(z) and B(z) are 
class invariants, and we can compute the class polynomials Hjj(A;x) and Hd(B;x) using 
standard algorithms such as those found in [5, 14, 16]. Under the GRH, we can compute these 
class polynomials in 0(|D|) expected time, which is negligible compared to the 0(|.D| 5 / 2 ) 
expected running time of Algorithm 2. 

For g = A,B, we can use H D (g; x) to uniquely determine a root gj, of ^ g (x, j k ) corresponding 
to a singular modulus j k by computing the unique root of the linear polynomial 

f k (g;x) := gcd(V g (x, j k ),H D (g;x)). 

This is useful because we would otherwise have 6 possible values of g k to choose from; in 
both cases ^ 9 {x,j k ) has degree 12, and 6 of its roots lie in the ring class field. In the context 
of Algorithm 2, we can use the values g k to replace the quantities r n (Jk) in step 3.c.ii that 
require A n (z) to be a rational function of j(z). For this purpose we let a k and b k denote 
the unique roots of the polynomials f k (A;x) and f k (B;x), respectively, and let a k denote 
a k /{3k{jk - 1728)). 

There is one other issue to consider. The coefficients of the class polynomials Hd(A;x) 
and He>(B;x) are not rational integers; they are algebraic integers in the quadratic field 
K = Q(v^D). This presents a potential difficulty for the CRT approach; while we always 
work modulo primes p for which D is a quadratic residue, we must make an arbitrary choice 
for the square root of D modulo p, and there is no clear way to make these choices consistently 
across many primes p. The following lemma implies that it does not matter which choice we 
make, we will get the same answer in either case. 

Lemma 3.8. Let P(z) = A(z)/(j(z)(j(z) — 1728) + B{z) r ){z) be as above. Then for any 
discriminant D = 1 — 24n, we have 

{P(a Q ) : a Q G Q»%/T (6)} = {Pfaj : a Q G Q^/r (6)}. 

Proof. For any modular function / as P (e.g. any weak Maass form) with real coefficients, we 
have 

{(/ | W 6 )(a Q ) : a Q G 0^/^(6)} = {/foj : a Q G Q^/r (6)}. 
Using the invariance of P under We, we get 

{P(a Q ) : a Q G Q^ m ,/r (6)} = (Ffaj : a Q G Q^/r (6)}. 

□ 

We now give the algorithm to compute the partition polynomial H^ &rt (x). 
Algorithm 3 

Input: A positive integer n. 

Output: The partition polynomial H^ avt (x) G Q[x]. 

1. Let 1 — 24n = v 2 D , where D is a fundamental discriminant. 

2. For each divisor u of v, let D = u 2 Dq and compute Hd(P; x) as follows: 
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a. Compute the class polynomials Hd(A;x) and Hr>(B;x). 

b. Using the height bound B P (D) defined below, compute the polynomial H D (P; x) using 
a modified version of Algorithm 2 in which r (jk) is replaced by = ak/(jk(jk ~ 1728)) 
and rx(jk) is replaced by bk, where and bk are as defined above. 

3. Compute H% art (x) G Q[x] via Lemma 3.7. 

The height bound Bp(D) is defined as 
(3.12) B P (D) := CiBjip) + h{D) log 

where Bj(D) is an explicit bound on the height of the Hilbert class polynomial derived as in 
Lemma 8 of [33]. Here C\ denotes an effectively computable positive constant. 

Remark. We have not tried to obtain the optimal constant for c\. However, it is reasonable to 
suspect that we can take c\ := 7/3, which we note is equal to degj(^^)/ deg x ( 1 i r ^) = 28/12 
(which dominates degj^s)/ deg_ x -(^ r s) = 18/12). 

Lemma 3.9. For all discriminants D = 1 mod 24 we have ht(H D (P;x)) < B P {D). 

Proof. The proof is analogous to the proof of Lemma 3.6, which in turn follows as in the 
proof of Lemma 8 of [33]. Decorating that proof with the asymptotic properties of the Fourier 
expansion of the function P(z), which is the image of a simple weight —2 weakly holomorphic 
modular form under the differential operator <9_o := tt^ ' 4- H t 1 ; ; i gives the desired 

V z 2m dz 27rlm(z)' to 

result. □ 

We now prove Theorem 1.3 given in the introduction, which we restate here. 

Theorem 1.3. For all positive integers n, Algorithm 3 computes H^ rt {x). Under the GRH, 
its expected running time is 0(n 5 / 2 ) and it uses 0(n 2 ) space. 

Proof. The correctness of Algorithm 3 follows from Lemmas 3.7 and 3.8, and the correctness 
of Algorithm 2. For the complexity bound, we first note that the degree of H^ avt (x) is given 
by the Hurwitz-Kronecker class number H(l — 24n) = h(D), where D varies over negative 
discriminants dividing 1 — 24n. It is known that H(D) = 0(h(D) (log log |-D|) 2 ); see, e.g., [33, 
Lemma 9]. The complexity bounds then follow from the height bound in Lemma 3.9 and the 
complexity bounds in Theorem 1.2. □ 

To simplify the practical implementation of Algorithm 3, we make the following remark: it 
is not actually necessary to compute the class polynomials Hd(A; x) and He>(B; x). Instead, 
for each singular modulus jk, one can simply compute all 36 possible combinations + tijk, 
where si, . . . , s 6 are the roots of ^^(x,j^) that lie in ¥ p (where p = 11 mod 12 is a suitable 
prime) and t±, . . . ,te are the roots of ^b{x, jk) that lie in ¥ p . For all but finitely many primes 
p, exactly 32 of these 36 values will be distinct, and there will be two pairs of repeated values, 
corresponding to Pk = a k + b k jk and — P k = —a k — bkjk- We do not prove this claim here, but 
observe that Lemma 3.8 guarantees that the value Pk will be repeated, so if one in fact finds 
the situation modulo p to be as claimed (exactly two pairs of repeated values that differ only 
in sign), then the end result will be provably correct. For the handful of primes p where the 
claim does not hold, one simply discards p and selects another suitable prime in its place. 
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Using the observation above, one may compute the polynomial 

f{x) = \D\ 2h ^ H(x 2 -P%) = (-l) h ^\D\ 2h ^H D (P;x)H D (P ] —x) mod p 

k 

modulo a sufficient number of primes p (using a suitably increased height bound), and then 
apply the CRT to obtain the integer polynomial f(x) which may then be factored in 1i[x] to 
yield the required polynomial Hd(P;x). 

4. Numerical examples 

As a first example, let us compute iff ar (x) using Algorithm 3, recapitulating the example 
given in the introduction of [8]. We have n = 1, and the discriminant D = 1 — 24n = —23 is 
fundamental, sou= 1. We begin by computing the class polynomials 

H-2s{A;x) = x 3 + (264101659831625A- 76898070951625)/2x 2 

+ (4866595720359935196250A+ 1237728700002625503750) x 

+ (-14048754886813637262794029921875A - 31056014444792221417574181765625) /2, 
H_ 23 (B;x) = x 3 + (-35487375 A - 35487375) x 2 + (6837889760625 A - 75216787366875)/2 x 
+ (842331597312734375A+ 2863927430863296875)/2, 

where A denotes a square root of —23. We then compute the heuristic height bound 

B P (-23) = 7 /a Bj(-23) + h(-23) log 23 « 83.25 

for H-2z{P K ,x), using [33, Lemma 8] to compute Bj(— 23) ~ 31.65. 

Remark. We refer to this bound as a heuristic bound because we used Ci = 7/3 as in the 
discussion in the previous section. Implementing the algorithm with this bound for n < 750 
always gave the correct values for p(n). Moreover, in every case the polynomials Hd(P;x) 
computed by the algorithm split into linear factors over the ring class field for the order for 
discriminant D. 

We now use Algorithm 2 to compute iJ_ 23 (P;a;), with r (j k ) = a k /(j k (j k — 1728)) and 
r i(jk) — b k . The order O of index 11 in the order of discriminant —23 is suitable for the 
integer 23, and the primes in the set 

S = {1562207, 2744591, 4294607, 6454031, 7089107, 10010291} 

are all suitable for O and 23, with ripes^SP ~ 91.93 > Bp(— 23) + log2. We then compute 

#_ 23 (j ; x) = x 3 + 3491750x 2 - 5151296875s + 12771880859375 

using [33, Alg. 2]. 

Starting with the first prime p = 1562207, we compute $23 mod p using Algorithm 1.1 
(which in this case just calls [6, Alg. 2.1], since 23 is prime). We then find the roots j k of 
-^-23(2) mod p, and for each j k we compute: 

• 7* = (2/3 ,2 - #i,i)/A>,i mod V (via (2.9), using $23 mod p and j k ). 

• a k as the unique root of f k (A; x) = gcd(f/_ 2 3(A; x), ^^(x,j k )) mod p. 

• b k as the unique root of f k (B;x) = gcd(if_ 2 3(-B; x), ^b{x, jk)) mod p. 

• P k = a k / (jfc(1728 - j k )) + b k ^ k . 
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For the prime p = 1562207 the results of these computations are summarized below. 





k = 1 


k — 2 


k — 3 


jk 


244476 


467416 


482979 


Ik 


1461486 


587848 


220836 


a k 


1201792 


98544 


239915 


b k 


1120135 


560362 


531933 


P k 


1352290 


519913 


1252234 



Using the constants C\ = c 2 = 1 in our definition of a good modular function, we compute 
f(x) = 23 3 Y[(x - P k ) = 12167s 3 + 1282366a; 2 + 337961s + 1150855 (mod 1562207) 

k 

as the reduction of \D\ h ^Hc(P;x) modulo p. Repeating this process for the remaining 
primes in S yields the following polynomials \D\ h ^ D 'Hu(P'x) mod p. 

12167a; 3 + 2464750a; 2 + 1900168a; + 391209 (mod 2744591) 

12167a; 3 + 4014766a; 2 + 1900168a; + 3491241 (mod 4294607) 

12167a; 3 + 6174190a; 2 + 1900168a; + 1356058 (mod 6454031) 

12167a; 3 + 6809266a; 2 + 1900168a; + 1991134 (mod 7089107) 

12167a; 3 + 9730450a; 2 + 1900168a: + 4912318 (mod 1001029) 

Applying the Chinese remainder theorem, and using the fact that Yl P es P * s more than twice 
the absolute value of the largest coefficient of 23 3 if_23(-P; x), we obtain 

23 3 #_ 23 (P; x) = 12167a; 3 - 279841a; 2 + 1900168a; - 5097973 E Z[x]. 

We note that the coefficients of the above polynomial are all much smaller than the product 

Y[p = 5398465666938830659283417535896039, 

pes 

indicating that our height bound Bp(— 23) is actually bigger than it needs to be. Dividing by 
23 3 and applying Lemma 3.7 we obtain 

3599 

Hr\x) = H^ 23 (P; x)=x 3 - 23a; 2 + -^x - 419, 

which completes the execution of Algorithm 3 for n = 1. As proven in [8], if we divide the 
trace of H^ art (x) by 24n — 1, we obtain the nth partition number p(n). In this case we have 
23/23 = 1 =p(l). 

We now consider the case n = 24, which is the least n for which 1 — 24n is not a fundamental 
discriminant. We have 1 — 24 • 24 = —575 = — 5 2 • 23, and Lemma 3.7 then implies that 

(4.1) Hlf{x) = -H^ 23 (P; -x)H„ 575 (P; x). 
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We have already computed H_2z{P; x), and in the same way we may compute 

ff_ 575 {P; x) = x 18 - 905648 x 17 + 7864919720287/23 x 16 - 62085428963462224 a; 15 

+ 2500819220800663290310031/529 x 14 - 145570369368132345878793951/23 a; 13 

+ 758005997309239141979280480729944052789478182183267952 /3700897225 x 
- 274989755819545226019386671943056995003866543720439419/18504486125. 

From (4.1) we obtain 

HH r \x) = x 21 - 905625 a; 20 + 341932201569 a; 19 - 62077564185180110 a; 18 

+ 2500063855637055742916679/529 a; 17 - 143069773154897117981992275/23 a; 16 

- 248682508073724592034185083695904/60835 x 15 
+ 4721274513295479628753048946698042/2645 a; 14 

- 684240866701755248448205419660018178147/1399205 a; 13 
+ 828297525091153912001188772487055395656/12167 a; 12 

- 32704304695374273471069347508729088366971453/6436343 a; 11 
+ 290553028842402057481729080665422874771306601/1399205 a; 10 

- 15618334996574598433984982031615985504271825288372 /3700897225 a; 9 
+ 2971138261271289839650966142959376571788416952712/160908575 a; 8 

+ 67822191247241980381807708488720865403444300542792174/85120636175 a; 7 

- 10287891953477631667871642653944942982233172929865507/740179445 a; 6 

+ 120072172960067820695115892912976299403813878193923504758/1957774632025 a; 5 
+ 9442155332145807613622010202526881668517330792046133529/17024127235 a; 4 

- 944566531689753532003676376487531915501990271825184156855477/225144082682875 a; 3 

- 512515146501467199140764542151150418963279118308213518346717/9788873160125 a; 2 

- 35536755777441881604409993038352893457607117583456947874072/425603180875 a; 
+ 115220707688389449702123015544140880906620081818864116561/18504486125. 

If we divide the trace 905625 of (x) by 24 • 24 — 1 = 575 we obtain the partition number 
p(24) = 1575, as expected. A table of all the polynomials H^ art (x) for n < 750 is available 
online at http://math.mit.edu/~drew. 
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